library(haven)
library(ggplot2)
library(tidyverse)
library(patchwork)
library(kableExtra)
library(stargazer)
library(lmtest)
library(plm)
library(caschrono)
library(astsa)
library(forecast)
library(DT)
library(knitr)
library(plotly)
library(dplyr)
data <- read_dta("data_V3.dta")
data <- data[data$age != 996 & data$educyr <= 90,]

EXERCICE 1 : Les dépenses en soins dentaires

Le but de cet exercice est d’étudier les déterminants des dépenses en soins dentaires. Pour ce faire, nous disposons d’une base de données de panel. Elle a été construite à partir des réponses de 73 017 individus vivant aux Etats-Unis et qui ont été interrogés en 2015 et 2019.

Les variables sont les suivantes :

  • year : l’année d’enquête.
  • id : l’identifiant individuel.
  • dvexptot : le total des dépenses en soin dentaire.
  • age : l’âge, en années.
  • sex : le genre (2 pour une femme, 1 pour un homme).
  • educyr : le nombre d’années d’éducation.
  • inctot : le revenu total personnel.
  • hinotcov : n’a pas d’assurance santé (en excluant les plans d’épargne personnels) (1 non, 2 oui).
  • workev : a déjà travaillé (1 non, 2 oui, 7 a refusé de répondre, 8 non vérifié, 9 ne sait pas).
  • himcare : a l’assurance medicare (1 non, 2 oui).

Avant de commencer notre travail, nous avons étudié longuement la base de données afin de supprimer les variables aberrantes telles que des individus qui avaient plusieurs centaines d’années (ou bien certains qui avaient changé de sexe). Ce sont, a priori, des erreurs de saisie. Nous avons également modifié certaines variables afin de faciliter notre étude.

Question 1: Introduction

En quoi est-ce intéressant d’analyser les déterminants des dépenses en soin dentaire ?

Il est possible d’analyser les déterminants des dépenses en soins dentaires pour déterminer s’il existe des relations avec divers facteurs tels que le revenu, le sexe, le niveau d’éducation, l’âge, la possession d’une assurance maladie et la présence d’une assurance Medicare. Cette analyse peut fournir des informations sur les inégalités dans l’accès aux soins dentaires. Elle peut également aider à comprendre comment les différences socio-économiques influencent les dépenses en soins dentaires.

Question 2: Statistiques descriptives

Faites un tableau de statistiques descriptives avec les variables dont vous disposez.

Les individus dans notre base de données ont en moyenne 48 ans et sont tous des adultes. Ils ont généralement passé entre 12 et 16 années à étudier. Leurs revenus sont très variés, ce qui nous indique que les individus proviennent de plusieurs classes sociales différentes. En moyenne, les dépenses en soins dentaires sont de 307,5 dollars par an et par personne, mais cette moyenne est influencée par des valeurs extrêmement élevées. La médiane de ces dépenses est de 0 dollars, ce qui semble cohérent car peu de personne ont des dépenses dentaires, mais lorsqu’une personne doit payer pour des soins dentaires, cela coûte cher.

summary(data[,c("age","educyr","inctot","dvexptot")]) %>%
  kbl() %>%
  kable_minimal()
age educyr inctot dvexptot
Min. :19.00 Min. : 0.00 Min. :-309948 Min. : 0.0
1st Qu.:33.00 1st Qu.:12.00 1st Qu.: 10800 1st Qu.: 0.0
Median :48.00 Median :13.00 Median : 25021 Median : 0.0
Mean :48.38 Mean :13.02 Mean : 36915 Mean : 307.5
3rd Qu.:62.00 3rd Qu.:16.00 3rd Qu.: 50000 3rd Qu.: 200.0
Max. :85.00 Max. :17.00 Max. : 731653 Max. :81000.0

Il y a 54% de femmes et 46% d’hommes dans la base de données. Il y a donc un peu plus de femmes que d’hommes, mais le déséquilibre n’est pas très important. Selon les données de la Banque Mondiale, la part de femme aux USA varie entre 50 et 51% entre 1960 et 2021. (source).

Un quart des individus présent dans la base de données ont l’assurance medicare, c’est un programme de santé publique visant les personnes handicapées et les personnes âgées afin de leur donner une couverture medicale.

Nous avons 11% des individus qui n’ont pas d’assurance santé.

(ggplot(data, aes(x = sex)) +
  geom_bar(color = "black", fill = "#267334", alpha = 0.8, width=0.6) +
  ggtitle("Sexe") +
  labs(subtitle = "1 : homme et 2 : femme")+
  ylab("Effectifs")+
  theme_bw()+ 
ggplot(data, aes(x = year)) +
  geom_bar(color = "black", fill = "#267334", alpha = 0.8) +
  ggtitle("Année") +
  ylab("Effectifs")+
  theme_bw())/
(ggplot(data, aes(x = himcare)) +
  geom_bar(color = "black", fill = "#267334", alpha = 0.8, width=0.6) +
  ggtitle("Assurance medicare") +
  labs(subtitle = "1 non et 2 oui")+
  ylab("Effectifs")+
  theme_bw()+ 
ggplot(data, aes(x = hinotcov)) +
  geom_bar(color = "black", fill = "#267334", alpha = 0.8, width=0.6) +
  ggtitle("N’a pas d’assurance santé") +
  labs(subtitle = "1 non et 2 oui")+
  ylab("Effectifs")+
  theme_bw()) 

Nous pouvons constater que 89% des individus présents dans la base de données ont déjà travaillé. On a aussi la surprise de voir qu’une petite poignée d’individus ne savent pas s’ils ont déjà travailler ou non.

tab <- as.matrix(table(data[,c("workev")]))
rownames(tab) <- c("non","oui","a refusé de répondre","non vérifié","ne sait pas")
tab <- cbind(tab,tab/sum(tab)*100) 
colnames(tab) <- c("Effectifs","Proportions (%)")
tab %>% round(2)%>%
  kbl(caption = "A déjà travaillé ?", booktabs = T, linesep = "") %>%
  kable_minimal(full_width = F)
A déjà travaillé ?
Effectifs Proportions (%)
non 11018 10.60
oui 92102 88.58
a refusé de répondre 568 0.55
non vérifié 25 0.02
ne sait pas 266 0.26

Par la suite, nous allons modifier la varible workev, nous garderons seulement les réponses oui et non, elle fonctionnera donc comme les autres indicatrices.

Question 3: Modèle poolé

Estimez le modèle poolé suivant en corrigeant les écarts-types

\(dvexptot_{it}= \beta_0 + \beta_1age_{it} + \beta_2sex_{it} + \beta_3educyr_{i} + \beta_4inctot_{it} + \beta_5hinotcov_{it} + \beta_6workev_i + \beta_7himcare_{it} + \upsilon_{it}\)

avec \(\upsilon_{it} = \epsilon_{it} + c_i\)

data$sex <- as.factor(data$sex)
data$himcare <- as.factor(data$himcare)
data$hinotcov <- as.factor(data$hinotcov)
data$workev <- as.factor(data$workev)
library(stats)

data <- data[data$age != 996 & data$educyr <= 90 & (data$workev == 1 | data$workev == 2),]
pool <- plm(dvexptot ~ age + sex + educyr + inctot + hinotcov + workev + himcare, data = data,
                   model="pooling", index=c("id","year"))

pool_cor <- coeftest(pool, vcov=vcovHC(pool, method="arellano"))

En corrigeant les écarts-types, nous avons amélioré la précision de nos estimations, à l’exception de la variable himcare.

De plus, nous constatons que presque tous les coefficients du modèle sont significatifs au seuil de 1%, ce qui indique que chaque variable explicative a un impact statistiquement significatif sur la variable dépendante. workev est la seule variable non significative.


Variable dépendante : dvexptot
Modèle poolé Modèle poolé corrigé
age 4.159*** 4.159***
(0.283) (0.283)
sex2 63.037*** 63.037***
(7.089) (7.492)
educyr 26.216*** 26.216***
(1.280) (1.264)
inctot 0.002*** 0.002***
(0.0001) (0.0001)
hinotcov2 -97.688*** -97.688***
(11.664) (7.426)
workev2 6.474 6.474
(11.752) (11.990)
himcare2 72.639*** 72.639***
(11.882) (14.296)
Constant -349.614*** -349.614***
(22.216) (22.070)
Observations 103,120
R2 0.024
R2 ajusté 0.024
F Statistic 363.255*** (df = 7; 103112)
Note: p<0.1; p<0.05; p<0.01

Notre individu de référence est un homme qui n’a pas d’assurance Medicare, qui a une assurance santé privée et qui n’a jamais travaillé. Le coefficient associé à himcare est positif, montrant que l’assurance Medicare peut augmenter les dépenses de santé dentaire. Le coefficient de hinotcov est négatif, reflétant l’impact négatif de ne pas avoir d’assurance santé sur les dépenses de santé. Les années d’études, l’âge, un revenu élevé et l’expérience professionnelle sont des facteurs qui augmentent les dépenses de santé. Etre une femme a également un impact positif sur les dépenses de santé.

Question 4: Modèle à effets aléatoires

Réécrivez l’équation 1 dans le cadre d’un modèle à erreurs composées (ou effets aléatoires). Estimez ce modèle (n’oubliez pas de corriger les erreurs-types).

\(dvexptot_{it}= \beta_0 + \beta_1age_{it} + \beta_2sex_{it} + \beta_3educyr_{i} + \beta_4inctot_{it} +\beta_5hinotcov_{it} + \beta_6workev_i + \beta_7himcare_{it} + c_i + \epsilon_{it}\)

pool_ea <- plm(dvexptot ~ age + sex + educyr + inctot + hinotcov + workev + himcare, data = data,
                   model="random", effect="individual",index=c("id","year"))
pool_ea_corr <- coeftest(pool_ea, vcov=vcovHC(pool_ea, method="arellano"))

Variable dépendante : dvexptot
Modèle EA Modèle EA corrigé
age 4.104*** 4.104***
(0.297) (0.284)
sex2 62.651*** 62.651***
(7.462) (7.476)
educyr 26.313*** 26.313***
(1.345) (1.262)
inctot 0.002*** 0.002***
(0.0001) (0.0001)
hinotcov2 -97.005*** -97.005***
(12.096) (7.553)
workev2 6.564 6.564
(12.315) (12.008)
himcare2 74.702*** 74.702***
(12.452) (14.360)
Constant -347.061*** -347.061***
(23.305) (22.091)
Observations 103,120
R2 0.022
R2 ajusté 0.022
F Statistic 2,284.897***
Note: p<0.1; p<0.05; p<0.01

Question 5: Comparaison Poolé & E.A.

Comparez les résultats à ceux obtenus par les MCO.


Modèle poolé Modèle EA
age 4.159*** 4.104***
(0.283) (0.284)
sex2 63.037*** 62.651***
(7.492) (7.476)
educyr 26.216*** 26.313***
(1.264) (1.262)
inctot 0.002*** 0.002***
(0.0001) (0.0001)
hinotcov2 -97.688*** -97.005***
(7.426) (7.553)
workev2 6.474 6.564
(11.990) (12.008)
himcare2 72.639*** 74.702***
(14.296) (14.360)
Constant -349.614*** -347.061***
(22.070) (22.091)
Note: p<0.1; p<0.05; p<0.01

Il y a peu de différence entre le modèle poolé et le modèle à effets aléatoires. Les signes, les écarts-types et les coefficients sont similaires dans les deux modèles. Cependant, le choix peut se porter sur un modèle ou sur l’autre en fonction des hypothèses de base. En effet, un modèle poolé ou un modèle à effets aléatoires doit respecter la condition \(Cov(c_i,x_it) = 0\). Mais on sait que l’estimateur des MCO sur le modèle poolé est convergent (si la condition précédente est vérifiée) mais non efficace, c’est pourquoi un modèle à effets aléatoires est préférable.

Notons que l’on peut considérer que le revenu total d’un individu peut être corrélé avec le nombre d’années d’éducation qu’il a reçues. Dans ce cas, on est un peu embêté pour choisir entre ces deux modèles car ils sont probablement non convergents. Cependant dans le cas d’un modèle à effets fixes nous pouvons avoir un modèle convergent dans le cas où \(Cov(c_i,x_it) ≠ 0\)

Question 6: Modèle à effets fixes

Ecrivez et estimez le modèle à effets fixes individuels. Comparez les résultats à ceux obtenus avec le modèle à effets aléatoires.

Nous allons effectuer un modèle within. Comme nous sommes dans un modèle à effets fixes on s’attend à ce que nous n’ayons pas de coefficients pour les variables qui ne varient pas dans le temps.

\(dvexptot_{it} - \overline{dvexptot_{it}} = \beta_1(age_{it} - \overline{age_{it}}) + \beta_2(sex_{it} - \overline{sex_{it}}) + \beta_3(educyr_{i} - \overline{educyr_{i}})+ \beta_4(inctot_{it} - \overline{inctot_{it}})\) \(+ \beta_5(hinotcov_{it}- \overline{hinotcov_{it}}) + \beta_6(workev_i - \overline{workev_it}) + \beta_7(himcare_{it} - \overline{workev_i}) + (\epsilon_{it} - \overline{\epsilon_{it}})\)

fixed <- plm(dvexptot ~ age + sex + educyr + inctot + hinotcov + workev + himcare, data = data,
                   model="within", index = c("id", "year"))

fixed_cor <- coeftest(fixed, vcov=vcovHC(fixed, method="arellano"))

Variable dépendante : dvexptot
Modèle EF Modèle EA
age -6.064 4.104***
(7.195) (0.284)
sex2 -39.470** 62.651***
(15.847) (7.476)
educyr 26.313***
(1.262)
inctot 0.001*** 0.002***
(0.0003) (0.0001)
hinotcov2 -37.540 -97.005***
(27.801) (7.553)
workev2 -19.919 6.564
(33.522) (12.008)
himcare2 242.199*** 74.702***
(80.002) (14.360)
Constant -347.061***
(22.091)
Note: p<0.1; p<0.05; p<0.01

Pour toutes les variables fixes dans le temps, nous n’avons donc pas de coefficient dans le modèle à effets fixes, contrairement à un modèle à effets aléatoires. De plus, nous constatons que dans le modèle à effets aléatoires, peu de coefficients sont statistiquement significatifs, seuls inctot et himcare sont significatifs au seuil de 1%, alors que dans notre modèle à effets aléatoires, presque tous les coefficients sont significatifs. Cependant, prenons la variable himcare, qui est significative dans les deux modèles, on constate qu’elle a un écart-type très important dans le modèle à effets fixes (le coefficient est également bien plus élevé).

Question 7: Comparaison EA & EF

Entre le modèle à effets aléatoires, et le modèle à effets fixes, quel modèle choisiriez-vous ? Expliquez la procédure que vous utilisez pour déterminer votre choix.


EA EF
age 4.104*** -6.064
(0.284) (7.195)
sex2 62.651*** -39.470**
(7.476) (15.847)
educyr 26.313***
(1.262)
inctot 0.002*** 0.001***
(0.0001) (0.0003)
hinotcov2 -97.005*** -37.540
(7.553) (27.801)
workev2 6.564 -19.919
(12.008) (33.522)
himcare2 74.702*** 242.199***
(14.360) (80.002)
Constant -347.061***
(22.091)
Note: p<0.1; p<0.05; p<0.01

Il n’y a donc apriori aucun doute sur le modèle à choisir, c’est le modèle à effets aléatoires car les coefficients sont significatifs, les écarts-types sont plus faibles et, en plus, nous avons des coefficients pour les variables fixes dans le temps. Cependant il ne faut pas oublier que pour un modèle à effets aléatoires, nous devons respecter cette hypothése \(Cov(c_i,x_it) = 0\) si nous voulons un modèle convergent. L’avantage d’un modèle à effets fixes c’est qu’il est convergent si \(Cov(c_i,x_it) ≠ 0\), et ici on peut penser que cette corrélation n’est pas nulle (corrélaltion entre l’éducation et le revenu notamment). Afin d’être sûr, nous pouvons réaliser un test d’Hausman, cela nous permettra de nous aider à faire notre choix entre le modèle à effets fixes ou à effets aléatoires.

\[\begin{cases} H_{0}: Cov(c_i,x_{it}) = 0\\ H_{1}: Cov(c_i,x_{it}) ≠ 0 \end{cases}\]

library(lmtest)
test<-phtest(pool_ea,fixed)
tabtest <- rbind(c("Statistique","$p-value$","$ddl$","Alternative"),
c(round(test$statistic[[1]],2), round(test$p.value,5),
test$parameter[[1]], test$alternative))

tabtest %>%
  kbl(caption = "Test d'Hausman", booktabs = T, linesep = "", ) %>%
  kable_minimal(full_width = T)
Test d’Hausman
Statistique \(p-value\) \(ddl\) Alternative
28.77 7e-05 6 one model is inconsistent

On rejette donc \(H_O\) car la p-value est inférieur à 5%. Donc selon le test d’Hausman \(Cov(c_i,x_it) ≠ 0\). On va donc préferer un modèle à effets fixes afin d’avoir un modèle convergent.

Question 8: Le meilleur modèle

Parmi maintenant l’ensemble des modèles présentés, quel est le modèle que vous choisiriez ?

On va alors préférer choisir le modèle à effets fixes. En effet, selon le test réaliser précedemment il y aurait une corrélation entre les effets individuels et nos variables observés.


Dependent variable:
poolé EA EF
age 4.159*** 4.104*** -6.064
(0.283) (0.284) (7.195)
sex2 63.037*** 62.651*** -39.470**
(7.492) (7.476) (15.847)
educyr 26.216*** 26.313***
(1.264) (1.262)
inctot 0.002*** 0.002*** 0.001***
(0.0001) (0.0001) (0.0003)
hinotcov2 -97.688*** -97.005*** -37.540
(7.426) (7.553) (27.801)
workev2 6.474 6.564 -19.919
(11.990) (12.008) (33.522)
himcare2 72.639*** 74.702*** 242.199***
(14.296) (14.360) (80.002)
Constant -349.614*** -347.061***
(22.070) (22.091)
Note: p<0.1; p<0.05; p<0.01

Question 9: Conclusion

Comment synthétiseriez vous vos résultats concernant les dépenses en soin dentaire ?

Avec la décision prise à la question précédemment on peut dire que les dépenses en soin dentaire s’explique par trois variables, le sex, la possesion ou non de l’aide Medicaire et du revenu des individus. Plus un individu est riche et plus il pourra allouer d’argent à ses dépenses dentaires. L’aide medicare permet également d’augmenter les dépenses en soin dentaire.

Cependant on peut-être interloqué de la présence du sexe comme variable explicative alors que nous avons décider de prendre un modèle à effets fixes.

Question 10: Données incohérentes

Il est possible que les réponses croisées de certaines variables soient incohérentes pour certains individus. Le cas échéant, repérez ces observations et, si nécessaire, estimez les différents modèles en les incluant puis en les supprimant. Vos résultats changent-ils ?

On remarque que dans le modèle à effets fixes, il y a encore la variable sex, cela signifie alors que certains individus changent de sexe au cours du temps. Cela semble venir d’une erreur de remplissage de la base de données, nous pouvons retirer ces individus et observer s’il y a un changement sur le modèle.

progressite <- data %>%
  group_by(id) %>%
  summarize(n_unique_sex = n_distinct(sex)) %>% filter(n_unique_sex > 1) 

data$SexC <- ifelse(data$id == 21609 | data$id == 45492 | data$id == 48018, "Change","NoChange")
data <- data[data$SexC == "NoChange",]
data <- data[,-11]
fixed2 <- plm(dvexptot ~ age + sex + educyr + inctot + hinotcov + workev + himcare, data = data,
                   model="within", index = c("id", "year"))

fixed_cor2 <- coeftest(fixed2, vcov=vcovHC(fixed, method="arellano"))

Variable dépendante : dvexptot
Modèle EF Modèle EF corrigé
age -6.063 -6.063
(7.167) (7.195)
inctot 0.001*** 0.001***
(0.0002) (0.0003)
hinotcov2 -37.533 -37.533
(33.388) (27.801)
workev2 -19.920 -19.920
(71.695) (33.522)
himcare2 242.199*** 242.199***
(62.543) (80.002)
Observations 103,114
R2 0.001
R2 ajusté -1.369
F Statistic 6.435*** (df = 5; 43502)
Note: p<0.1; p<0.05; p<0.01

Hormis le fait que le coefficient associé à la variable sex n’est plus présent, nous n’observons pas de changement particulier que ce soit sur les signes des coefficients, leur valeur ou bien même les écart-types. Désormais les variables qui expliquent les dépenses dentaires sont inctot & himcare.

EXERCICE 2 : Théorie

Question 1: Modèle ARIMA

Simulez le modèle suivant

\[Y_t=0.5Y_{t-1}+0.6\epsilon_{t-1}+\epsilon_t \text{ avec } T = 200\]

set.seed(1)

AR1 <-arima.sim(list(ar= 0.5, ma=0.6), n=200)

ggplotly(ggplot(data.frame(AR1), aes(x = 1:200, y = AR1)) +
  geom_line(color = "#577559") +
  ggtitle("Modèle ARIMA(1,0,1)") +
  xlab("Temps") +
  ylab("Valeur") +
  theme_bw()) 

Question 2: Fonction ACF et PACF

Représentez la fonction d’autocorrélation, et la fonction d’autocorrélation partielle. Commentez-les.

acf2y(AR1, 10) %>% 
  kbl( booktabs = T, linesep = "") %>%
  kable_minimal(full_width = T)
LAG ACF1 PACF
1 0.7155812 0.7155812
2 0.3465992 -0.3390911
3 0.1800755 0.2098042
4 0.1199695 -0.0685944
5 0.1038531 0.0838758
6 0.0813634 -0.0506896
7 0.0536202 0.0294945
8 0.0319742 -0.0203796
9 0.0482315 0.0833410
10 -0.0022429 -0.1973629

Pour l’ACF, on peut observer un pic à 0.7 environ puis une décroissance, cela veut dire qu’on a une corrélation significative entre les observations, donc le temps présent est corrélé à la dernière période passée ce qui peut indiquer la présence d’un modèle AR dans notre série chronologique.

Pour l’autocorrélation partielle, on a un pic positif à 0.7 ce qui peut indiquer qu’il y a une forte corrélation entre le temps \(t\) et \(t-1\), et donc un modèle autorégressif AR serait adapté. On observe aussi un pic négatif significatif à 0.3 environ, on peut donc penser qu’il y a une relation négative entre la période \(t\) et la période \(t-2\), ce qui indiquerait la présence d’un modèle à moyenne mobile (MA).

On peut donc conclure qu’il y a probablement une combinaison de modèles AR et MA dans notre série temporelle. Les pics positifs significatifs dans l’ACF indiquent une forte corrélation entre les observations passées et les observations actuelles, ce qui suggère un modèle AR. De même, le pic négatif significatif dans la PACF suggère une relation négative entre les observations actuelles et les observations passées, ce qui mène à un modèle MA. Par conséquent, le modèle ARMA pourrait être un bon choix pour modéliser cette série temporelle.

Question 3: modèle ARMA, AR et MA

Ajustez un ARMA(1,1), un AR(1) et un MA(1).

Modèle ARMA(1,1):

\[Y_t=0.18 +0.47Y_{t-1}+0.61\epsilon_{t-1}+\epsilon_t\]

ARMA <- arima(AR1,order = c(1,0,1))
ARMA
ar1 0.469***
(0.080)
ma1 0.608***
(0.083)
intercept 0.178
(0.203)
Observations 200
Log Likelihood -275.272
sigma2 0.913
Akaike Inf. Crit. 558.544
Note: p<0.1; p<0.05; p<0.01

Modèle AR(1):

\[Y_t=0.20 + 0.72Y_{t-1}+\epsilon_t\]

AR <- arima(AR1,order = c(1,0,0))

AR
ar1 0.719***
(0.049)
intercept 0.197
(0.261)
Observations 200
Log Likelihood -293.698
sigma2 1.100
Akaike Inf. Crit. 593.396
Note: p<0.1; p<0.05; p<0.01

Modèle MA(1):

\[Y_t=0.17 + 0.88\epsilon_{t-1}+\epsilon_t\]

MA <- arima(AR1,order = c(0,0,1))

MA
ma1 0.882***
(0.039)
intercept 0.166
(0.135)
Observations 200
Log Likelihood -287.852
sigma2 1.034
Akaike Inf. Crit. 581.704
Note: p<0.1; p<0.05; p<0.01

Question 4: Bruits blancs ?

Les résidus sont-ils des bruits blancs pour les trois modèles ?

On va d’abord effectuer un test de Box-Pierce pour vérifier si les résidus sont des bruits blancs donc qu’ils soient indépendants et pas corrélés.

\[\begin{cases} H_{0}: \epsilon_t \text{ sont indépendants (pas d'autocorrélation dans les résidus)}\\ H_{1}: \epsilon_t \text{ ne sont pas indépendants (autocorrélation dans les résidus)} \end{cases}\]

Si \(H_0\) est rejetée (p-value < \(\alpha\)), les résidus ne sont pas considérés comme des bruits blancs. Si on accepte \(H_0\) (p-value \(\geq \alpha\)), alors les résidus sont des bruits blancs.

Res <- ARMA$residuals
arm<-Box.test(Res)
ResAR <- AR$residuals
a<-Box.test(ResAR)
ResMA <- MA$residuals
m<-Box.test(ResMA)
tab<-rbind(c("ARMA", round(arm$p.value, 4)),
           c("AR", round(a$p.value, 4)),
           c("MA", round(m$p.value,4)))
colnames(tab) <- c("Modèle", "p-value")


tab %>%
  kbl(caption = "Comparaison", booktabs = T, linesep = "") %>%
  kable_minimal(full_width = T)
Comparaison
Modèle p-value
ARMA 0.8616
AR 5e-04
MA 8e-04

Les résidus du modèle ARMA(1,1) ont des p-value \(> 0.05\) donc c’est un bruit blanc. Etant donné que les p-value des modèles AR et MA sont inférieurs à \(0.05\), on en conclut que leurs résidus ne sont pas des bruits blancs.

On vérifie que les résidus proviennent d’une distribution normale avec un test de normalité (ici le Shapiro test). Si les résidus suivent une distribution normale alors les erreurs dans les modèles sont réparties de manière aléatoire et ne sont influencées par aucune tendance. Cela assure que les conclusions basées sur le modèle sont fiables et robustes.

\[\begin{cases} H_{0}: \text{les résidus suivent une distribution normale}\\ H_{1}: \text{les résidus ne suivent pas une distribution normale} \end{cases}\]

Si la p-value \(\geq 0.05\), alors on accepte \(H_0\) et les résidus suivent une distribution normale.

arm_shapiro <- shapiro.test(Res)
ar_shapiro <- shapiro.test(ResAR)
ma_shapiro <- shapiro.test(ResMA)

tab2 <- rbind(c("ARMA", round(arm_shapiro$p.value, 4)),
                  c("AR", round(ar_shapiro$p.value, 4)),
                  c("MA", round(ma_shapiro$p.value, 4)))
colnames(tab2) <- c("Modèle", "p-value")


tab2 %>%
  kbl(caption = "Comparaison", booktabs = T, linesep = "") %>%
  kable_minimal(full_width = T)
Comparaison
Modèle p-value
ARMA 0.4267
AR 0.8647
MA 0.1195

Ici, tous les modèles ont des résidus provenant d’une distribution normale.

Question 5: AIC

Quel autre critère pourriez-vous utiliser pour choisir entre les 3 modèles ? Expliquer son principe. Quel modèle choisiriez-vous entre les 3 modèles ?

L’autre critère qu’on pourrait utiliser est le critère d’informations d’Akaike (AIC).

Ce critère mesure la qualité de prédiction d’un modèle en comparant son erreur de prédiction aux informations apportées par son nombre de paramètres.

\[AIC = ln(\frac{SCR_{\epsilon}}{T}) + \frac{2(p+q)}{T}\]

tab3 <- rbind(c("ARMA", round(AIC(ARMA),2)),
           c("AR", round(AIC(AR),2)),
           c("MA", round(AIC(MA),2)))

colnames(tab3) <- c("Modèle", "AIC")

tab3 %>%
  kbl(caption = "Comparaison", booktabs = T, linesep = "") %>%
  kable_minimal(full_width = T)
Comparaison
Modèle AIC
ARMA 558.54
AR 593.4
MA 581.7

Plus l’AIC est faible et plus le modèle est de bonne qualité.

Le modèle qui a l’AIC le plus faible est considéré comme étant plus proche de la réalité, ici c’est le modèle ARMA.

Question 6: Cohérence

Est-ce cohérent avec le modèle de la question 1) ?

Le modèle choisit est cohérent car à la question 1 on utilise aussi un modèle ARMA qui est un modèle à la fois autorégressif et à moyenne mobile. Le fait que l’AIC montre que le meilleur modèle soit l’ARMA(1,1) vient renforcer la pertinence de ce modèle, car il prend aussi en compte les relations dans les données passées et les erreurs.

Le modèle ARMA est le plus adapté pour représenter notre simulation réalisée à la question 1, ce qui est logique étant donné que la simulation a été effectuée à partir d’un modèle ARMA.